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We search for translationally invariant states of qubits on a ring that maximize the nearest 
neighbor entanglement. This problem was initially studied by O'Connor and Wootters [Phys. Rev. 
A 63, 052302 (2001)]. We first map the problem to the search for the ground state of a spin 
1/2 Heisenberg XXZ model. Using the exact Bethe ansatz solution in the limit N — » oo, we 
prové the correctness of the assumption of O'Connor and Wootters that the state of maximal 
entanglement does not have any pair of neighboring spins "down" (or, alternatively spins "up"). 
For sufnciently small fixed magnetization, however, the assumption does not hold: we identify the 
rcgion of magnetizations for which the states that maximize the nearest neighbor entanglement 
necessarily contain pairs of neighboring spins "down" . 

PACS numbers: 03.67.Mn, 75.10.Pq, 03.67.-a 



I. INTRODUCTION 

The investigation of the role of entanglement in quan- 
tum and classical phase transitions, and morc gcncrally 
the role of entanglement in many-body quantum systems 
is one of the hottest interdisciplinary areas on the borders 
between quantum information, quantum òptics, atòmic, 
molecular, and condensed matter physics. Initially the 
studies of entanglement in many-body systems have been 
motivated by the possibility of employing entanglement 
for quantum computation in optical lattices 0, or pre- 
cision measurements with Bose-Einstein condensates . 
Recently, several lines of research have been followed: 

• Studies of local entanglement in spin systems 0, 
0, 0, 0- 0]- and more generally in various many- 
body systems (such as linear chains, see for instance 
[3. Isí ITolp. with particular attention to the role of 
entanglement in phase transitions. 

• Studies of the entropy of blocks of spins and the re- 
lated "area law" [Hj , indicating weak entanglement 
of blocks, and effects of criticality [H El El El- 

• Studies of localisable entanglement and entangle- 
ment correlation length that diverges at the criti- 
cal point 0, El and majorizes the Standard cor- 
relation length (see also 0). In particular, it 
has been shown that the localisable entanglement 
(cf . [ÏtI lïsl El ) is bounded from above by the en- 
tanglement of assistance [2(| and from below by 
correlation functions. It follows directly from these 
bounds that one can define an entanglement cor- 
relation length that diverges in quantum critical 
systems. 

• Studies of multipartite entanglement in many-body 
systems pi). 



• Studies of dynamics and generation of entangle- 
ment in many-body systems 0, IÜEÍEIIsI In 
particular, implementations of the "one-way quan- 
tum computer" and short range teleportation j2|J 
of an unknown state has been proposed b y using 
the dynamics of spin systems in Refs. [13, lïï I2M 

eiimei 

• Studies of quantum information and entanglement 
theory inspired numerical codes to simulate quan- 
tum systems |3fij . 

A different approach to the study of entanglement in 
many-body systems has been proposed in two papers by 
Wootters and coworkers In these papers, instead 

of looking at a specific Hamiltonian, the authors asked 
the fundamental question "what is the maximal entan- 
glement between two neighboring sites of an entangled 
ring with translational invariance?" Here, an entangled 
ring is a chain of spins with periòdic boundary condi- 
tions. Due to the so-called "monogamy of entanglement" 
it is impossible for a site to be maximally entangled with 
both its neighbors: shared entanglement is always less 
than maximal |32l l33| . In Ref. [3l| the question of the 
upper limit for the nearest neighbor (NN) entanglement 
was simplified by introducing two additional restrictions 
on the allowed states: 

(i) The state of N spins 1/2 is an eigenstate of the 
z-component of the total spin (i.e. it has a fixed 
number of "down" spins p < N/2) |43| . 

(ii) Neighboring spins cannot both be "down" . 

Obviously, one can equally well study the same prob- 
lem in terms of spins "up", when N > p > N/2. Both 
restrictions are based on an educated guess for the op- 
timal states for the general problem. O'Connor and 
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Wootters (OW) solved the restricted optimization prob- 
lem by relating it to an effective Hamiltonian for the 
one-dimensional ferromagnetic XY model, and found the 
maximal nearest-neighbor concurrence (cf. Sec. [nj for 
given N and p to be 
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For given N and p, Eq. Q provides a lower bound for 
the problem without restriction (ii). It may, or may not 
happen that C can be increased by also allowing states 
where two neighboring spins are "down" . We have re- 
cently studied finite size rings and found that for a fixed 
p restriction (ii) tends to play a less important role as 
N is increased |34j: For p close to N/2 one can increase 
the concurrence significantly by dropping restriction (ii), 
but for p < 7V/3 OW's result is the optimal one. In fact, 
already in Ref. 31] it was shown that for all even N the 
ground state of a Heisenberg spin 1/2 antifcrromagnetic 
ring maximizes the concurrence among the zero magne- 
tization (p = N/2) states although it violates restriction 
(ii). 

By optimizing with respect to p one obtains a lower 
bound on the overall optimal concurrence, i.e. without 
any restrictions bcsides the translational invariance. In 
the limit N — > oo, the optimal number of spins "down" 
in Eq. (JIJ approaches p op t « 0.301 N. This leads to an 
asymptotic value of m 0.434. Although Ref. 

as well as our previous work |34j showed evidence for 
optimality of this number, whether it can be improved 
was, so far, an open problem. 

Wolf, Verstraete, and Cirac have in Ref. [35[ directly 
related OW's type of problems of looking for translation- 
ally invariant states that maximizc local entanglcment to 
the study of the ground state of a suitably defined "par- 
ent" Hamiltonian. In this paper we use this method and 
employ the known exact solution of the corresponding 
parent Hamiltonian to prové rigorously that: 

(A) In the limit N — > oo the translationally invariant 
state that maximizes the NN entanglement with- 
out any restriction coincides with the state found 
by OW at the optimal value of p ~ 0.3017V. This 
means that it is not a superposition of states with 
different p vàlues and does not contain simultane- 
ously neighboring spin "up" and neighboring spin 
"down" pairs. 

(B) For fixed p sumciently close to N/2, i.e. for suf- 
ficiently small magnetizations, assumption (ii) is 
not correct: the states that maximize the nearest 
ncighbor entanglement necessarily contain simul- 
taneously pairs of neighboring "up" and "down" 
spins. In the limit N — > oo we identify rigorously 
an interval of p/N for which this is the case and 
show strong numerical evidence that this interval 
is optimal. 



Our paper is organized as follows. In Section [l·l]we ap- 
ply the method of Ref. [3^ and derive the corresponding 
parent Hamiltonian for a N qubit ring. In Section IIIII 
we show the connection with the "classical" papers of 
Yang and Yang on the XXZ model. In Section IIVI we 
discuss briefly the regimes of parameters of interest and 
show that the present problem concerns the "difficult" 
parameter region of the phase diagram. In Section Ivlwe 
present the analysis based on the limit N — > oo of the 
Bethe ansatz solutions. We derive here the bàsic integral 
equation, the solution of which allows to calculate the 
desired energy of the system in question. In Section IVII 
the numerical results are discussed. In Section IVIII wc 
rigorously prové that the states that were conjectured 
in Ref. |3j] to maximize the NN entanglement and con- 
firmed by us, indeed provide the maximum of the NN 
entanglement for sufficiently small vàlues of p. We iden- 
tify the region of p/N where the latter statement does 
not hold. We conclude in Section IVIIII The short ap- 
pendix contains simple analytic bounds on optimal mag- 
nètic field for which the NN concurence is maximal. 



II. VARIATIONAL CONCURRENCE FORMULA 

In this paper we will use the concurrence as our en- 
tanglement measure. The concurrence [3(| is defined as 
C(p) = max{Ai — A2 — A3 — A4,0}, where Ai > A2 > 
A3 > A 4 > are the square roots of the eigenvalues of 
pp, and p = [cr y ® <Jy)p*(p y ® <jy) is the spin-flipped den- 
sity matrix. The optimization problems that we consider 
are complicated by the nonlinearity of the concurrence as 
function of the density matrix. In our previous work |34| , 
we showed how the optimization problem with fixed p can 
be reformulated as finding the ground state energy for 
each in a family of spin-chain Hamiltonians. This family 
is parameterized by a single real parameter and the opti- 
mal concurrence is minus the lowest ground state energy 
that oceurs when this parameter is varied. In this way 
a complicated non-linear problem in a high-dimensional 
space is replaced by a series of linear problems and one 
final one-parameter optimization. 

The derivation in Ref. [34| did not cover the case where 
superpositions of states with different p are allowed. To 
treat that case, we turn to Ref. |35j where the following 
general formula for the concurrence for systems of two 
spins 1/2 has been derived: 



C{p) = max jo, - d inf ^ tr [p (X ® X f F)] j 



(2) 



Here X is an arbitrary 2x2 matrix of determinant 1, 
while F is the flip (or swap) operator, interchanging the 
two qubits: 

f = ITTXTTI + II4XIII + IITXUI + IUXITI- (3) 

A useful parameterization of X is obtained using the sin- 
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guiar value decomposition |37| : 



X = U 



it 



(4) 



where t € [-00,00], U,V € U(2) and dett/ • detVt = 1. 
In fact, from Eq. J5J it is clear that we can restrict to 
U, V S ST/(2). We now rewrite 



X® X f F 

= (17 ® V) 
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Let us define the matrix in square brackets as h(t 2 ), i.e., 



h(s) = «ITTXTTI + -IUXUI - IITXTII - líDUTI- 

where s — t 2 . Then we can rewrite Eq. as 



(6) 



{U®V) f p(U®V)h(s)} \. (7) 



max { 0, — inf tr 

s,U,V 



Our goal is to maximize the concurrence over all p 
that can occur as nearest neighbor density matrices on a 
translationally invariant ring. If we always had U = V, 
it is easy to see that we could drop the infimum over 
U £ SU{2) in Eq. since if p = tr 3 ...Arr with T trans- 
lationally invariant then {U ®U)' S p{U ®U) = tr 3 ...jv(t^ <8 
. . .®Uyr(U®. . .®U) where (U®. . .®U)^T(U®.. .®U) 
is translationally invariant as well. To do the same for 
U =/= V, we can use the fact that h(s) is symmetric in the 
two qubits: 



tr Uu ® V) f p(U®V) h(s)] = tr [ph(s)] 



(8) 



where 



P = 



- { (U ® V)** p (U ® V) + (V ® U) 1 p (V ® U) } 



(9) 



If N is even and p = tr3 . jvT then p is a nearest neighbor 
density matrix belonging to the following translationally 



invariant state: 

{[U ®V]® ...®[U ®V]) ] x 

x T x ([U ®V}® . . . ® [U ®V]) 
+ ([V®U}®...®[V®U})ï x 

xT x {[V ®U}® ...®[V ®U}) 



(10) 



If N is odd, the above construction does not work: we 
cannot fit an integer number of U ® V terms on the ring. 
By placing as many terms U ® V as possible, and taking 
the translationally variant mixture of the resulting state, 
p can be approximated up to a factor 1/./V. In the limit 
of N — > 00 we can ignore this correction. 



III. THE PARENT HAMILTONIAN 

In this section we will follow the approach of Ref. 35] 
to derive the parent spin 1/2 XXZ Hamiltonian, i.e. the 
Hamiltonian whose ground state maximizes the NN con- 
currence. We also make the connection to the classical 
papers on the XXZ model by Yang and Yang [3^139114(1 ]. 

In the previous section we showed that in the limit 
N —> 00 

C max (Af) = maxC(p) = -mïtr[ph(s)} (11) 

P s,p 

where p = tr3...A/T for some translationally invariant T of 
N spins. The two-spin Hamiltonian can be written 
in terms of Pauli matrices as 

h(s) = S - (1 + a z ) (1 + a z ) + i (1 - a z ) (1 - a z ) 



— CT + (T + — <7_<7_ 

= \ { !J ' n · ~ n ·' CÍ ' 1 ' ~, ( S ' ~ ) { °' 1 1 

1 / 1 



(12) 



where a± = a x ±ia y . Instead of working with p and h(s) 
we can use the translational invariance and work with T 
and a Hamiltonian for the whole ring obtained by taking 
(|T2jl for each NN pair: 



Hwo\í(s) 



1 N 

— Y { - <tX +1 - aÍaÍ +1 - A{s)aíal +1 

on t-^i L x x y y \ t z z 

i=l V íà ) 



where 
A( S ) 



U S + -s 



2H(s)4-A(s)}, 



(14) 
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We have then reformulated the overall optimization prob- 
lem as 

r ax (iV) = maxCCp) = - inf tr[r#woif (a)l, (15) 
p s,r 

where p is restricted to arise from a translationally invari- 
ant state of N spins while the optimal T can automati- 
cally be chosen so since ífwoif is translationally invariant. 

An important observation can be made from Eq. 1)150. 
namely that as iïwoif commutes with the z component 
of the total spin, in the considered limit of ./V — > oo OW 
were right when they made assumption (i): The optimal 
state can indeed be chosen to have a definite number of 
spins "down" and thus does not contain superpositions 
of states with different p vàlues. Conversely, from our 
previous work [24j we know that Eq. 1)15(1 is also vàlid for 
any fixed p, i.e. we can write C max (N,p) on the left-hand 
side when making the appropriate restrictions on T. In 
summary, for fixed p the maximal concurrence is given 
by: 

r"(JV,p) = -ME os [H Wo}í (s),p], (16) 

s 

where i?Gs[ífwoif(s),p] is the "ground state" energy of 
Hwoif(s) in the manifold of states with p spins "down" . 
The overall maximal concurrence is given by furthcr op- 
timization over p or, equivalently, by using unrestricted 
ground state energies: 

C max (N) = maxC max (Y,;p) 

(17) 

= - mf E GS [H Woii (s)}. 

Let us now describe the connection with the work of 
Yangand Yang. In their seminal papers Yang and Yang 
(^ElMEU study this anisotropic Heisenberg XXZ Hamil- 
tonian (see e.g. ^ÏJ f°r morè recent work): 

ÍÍYang = ~\ ]T + + A^a^ 1 } , (18) 

i 

and they define / = Iíhiat^oo /at, where /at is half the 
energy per spin in the ground state with a given number 
p of spins "down" : 

f N (A,y) = ^E GS (H Yws ,p). (19) 

Here y is the average magnetization: 

S/=^<£4> = 1-|. (20) 

Since p is a conserved quantum number, one can include 
a magnètic field along z and only shift the energy of each 
eigenstate. The translation of the results of Yang and 
Yang to our optimization problems is therefore, 

£Gs[-ffwoif(s),p] = 

2f N (A(s),y)-H(s)y-~A(s), (21) 




Figure 1: Phase diagram for the XXZ chain in a magnètic 
field. For large negative A and small magnètic field 7i the 
ground state has perfect anti-ferromagnetic ordering (AF) 
with vanishing magnetization. When the magnètic field is 
increased beyond a certain critical value (lower bold line in 
the plot), a non- vanishing magnetization develops. The mag- 
netization increases with increasing field until finally all the 
spins are aligned. This aligned phase (A) is entered at an 
upper critical field indicated by the upper fat line in the plot. 
We are interested in the properties of the chain along the 
dashed line in the plot, i.e. we deal with the phase with non- 
vanishing, but also non-saturated magnetization. 

with s > 0. Note that A 2 - H 2 = 1. To find C max (N,p) 
we should minimize Eq. (|21[1 over s while keeping y fixed 
at the value corresponding to p [cf. Eq. (|20|l ]. To find the 
overall maximal concurrence C max (iV) we should further- 
more minimize over y. 



IV. THE PHASE-DIAGRAM OF THE XXZ 
MODEL 

We are considering the XXZ model in the limit TV — > 
oo. The second Yang and Yang paper [3íj deals with 
the properties of /(A, y) in exactly this limit. The third 
paper [40J contains information about the magnètic prop- 
erties, i.e. it is highly relevant when we also vary y in 
order to find the optimal fraction of spins "down" . In 
order to understand the regimes of parameters we are in- 
terested in and relate them to the known properties of 
the model, it is useful to look at the phase diagram of 
this model, displayed in Fig. ^ 

Let us first identify the region of the phase diagram 
which belongs to our parent Hamiltonian (|13fl . From 
Eq. (|14fl it is clear that as s varies from to oo, we move 
on a hyperbola in the A-7í plane: The s — case corre- 
sponds to (— oo, oo), whereas at s = 1 we are at the point 
of closest approach and cross the A axis in (—1,0), and 
as s — > oo we move back to infinity, but this time with 
negative magnètic field. Comparing this with Fig. ^ it is 
then easy and not surprising to see that the s-hyperbola 
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lies exactly in the "difficult" region of the phase diagram, 
i.e., the part with neither perfectly aligned spins nor per- 
fect anti-ferromagnetic order (between AF and A in Fig. 

®- 

Since a change of sign of the magnètic field will only 

interchange the role of spin "up" and spin "down" , we 
can ignore the negative TL branch and focus on s € [0, 1]. 
Then each point on the curve corresponds to exactly one 
A and we can thus parameterize the curve by A instead 
of s. The optimization is then done over A with the 
magnètic field always given by Tí = VA 2 — 1. 



V. THE INTEGRAL EQUATION 

The Bethe ansatz basically consists of the assumption 
that the wave function can be written as a sum of plane 
waves with a limited number of terms. If we are looking 
for a state with p spins "down" , only p wave numbers 
are needed. For the XXZ chain the first Yang and Yang 
article shows that this is indeed enough to produce the 
ground state wave function j^. For our purposes we 
should note that Yang and Yang give explicitly the equa- 
tions one needs to solve in order to find the ground state 
energy. In the limit of N — > oo, the number of wave num- 
bers naturally becomes infinite and the equation to find 
them becomes an integral equation for the wave number 
density. This equation mathematically has a form of a 
so called Fredholm equation of the second kind. After 
some reparameterization this equation attains the form 
(Eq. [7a] in Ref. [Ü): 



Rf a ) = - — 

da 2tt 



b ç\n 



(22) 



The unknown function here is R, which is the reparam- 
eterized density of wave numbers. The other functions 
depend parametrically on A, and in terms of the param- 
eter A = cosh _1 (— A), they are explicitly given by: 



dp sinh À 

da cosh A — cos a ' 

de _ sinh2A 

df3 cosh 2A — cos(a — (3) 



(23) 
(24) 



Let us point out the importance of the integration limit 
b in Eq. (|22|l : Whcn varying b, we get solutions corre- 
sponding to different vàlues of y. In fact, y is given by: 



tt(1 -y)= I R{a)da. 

-b 



(25) 



Note, however, that R also depends on ò, so the con- 
nection is not very obvious. In praxis (i.e. when doing 
numèrics) one solves Eq. Ij22(l for a range of the parame- 
ter b in order to find the result for the wanted vàlues of y. 
If one wants to optimize some quantity with respect to y, 



however, this can equally well be achieved by optimizing 
with respect to b. 

We are not primarily interested in R (which describes 
the state), but in /, which is the energy. It is given by: 



/(A,y) = -- 



A sinh A 



2tt 



-6 



R(a)^da. (26) 
da 



Again, / is written as a function of y, but in praxis the 
dependence is via b. 



VI. NUMERICAL SOLUTION OF THE 
INTEGRAL EQUATION 

A possible way to solve Eq. I|22|) is to turn the integral 
into a sum so that it becomes a matrix equation. This 
is called the Nystrom method The best way to dis- 
cretize an integral is not always equally spaced points; 
very often it is much more efficient to use a Gaussian 
Quadrature. This means that we evaluate the integrand 
at M points and make a weighted sum with weights 
{wk}- The points and the weights can be easily found in 
e.g. Mathematica. In this way, Eq. l)22[| becomes: 



Rk = 6= - y^wiKuRi 



where 



while 



R k = R(a k ) , £ fc 



sinh A 



cosh A — cos ak 



K M = 



1 



sinh 2À 



27T cosh 2 A — cos(afc — ai )' 



(27) 



(28) 



(29) 



It is clear that Eq. I|27|l is a matrix equation and that 
solving it cannot be harder than inverting 1 + K where 
Km = wiKki (no summation over l). 

The advantage of using Gaussian Quadrature is that 
one does not need too many points to get a very good es- 
timate of the integral for any sensible function. What ex- 
actly a "sensible function" is depends on the exact Gaus- 
sian Quadrature rule used. We use the simple Gauss- 
Legendre rule, assuming that R is well approximated by 
a polynomial on the interval [—6,6]. This is reasonable 
here because (|23ll and l|24l) are well-behaved for the vàl- 
ues of A we will consider. The final matrix equation can 
be solved very rapidly on a small size computer. A mod- 
erate value of M, however, means that our knowledge of 
R is restricted to a rather crude sampling; fortunately 
this is not a problem, since y and / are themselves inte- 
grals, and so can be evaluated with the full accuracy of 
Gaussian Quadrature. 

To give the reader an idea about the numèrics, we note 
that a simple Mathematica program will work very well 
with M < 30. To produce a plot /(A, y) versus y for A 
not too close to —1 it takes about one minute. To plot 
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Figure 2: We plot — Egs[Hwo\í]> that is, for each A the opti- 
mal y is chosen. It can be seen that the maximum is attained 
in the limit A — » — oo and that the limiting value coincides 
with OW's result indicated by the dashed line. 



the function of main interest, Eq. (|21Jl optimized over p 
(i.e. y, i.e. b) also only takes a few minutes. In Fig.[21we 
present the results of a Fortran program, which is (not 
surprisingly) much faster than the initial Mathematica 
code. The results indicate that OW's assumption (ii) 
was correct: When we plot BGs[-Hwoif] as function of A, 
we see that the optimal value of A is reached at — oo, 
and in this limit OW's result is recovered. We conclude 
that these simple numerical results indicate that the state 
that maximizes the NN concurrence without any restric- 
tions (i.e. optimized over p/N , i.e. y) coincides with the 
OW state fulfilling assumption (ii) (no NN pairs of spins 
"down"). 



VII. PERTURBATIVE CÀLCUL ATION 

Looking at Eq. I|21|) above we see that the finite value 
in Fig.|5]in the limit A — > — oo is obtained because some 
diverging terms happen to cancel each other. This is 
of course a great concern when doing numèrics since it 
means that a good relative precision (knowing the result 
to e.g. 1 ppm) may not be enough. The obvious strategy 
is to extract the solution in the strict limit A — > — oo. In 
this section, we will present a perturbative calculation in 
1/A. 

In zeroth order of the perturbation series we set 
coshA = oo in Eq. (|22|) . and arrive at the simple equa- 
tion: 



Ro(a) = l - J- / Ro(/3) 



d(3. 



(30) 



fmd the constant solution: 

Ro = - r and 



1 



yo = 



l+o l+o 

This means that in this limit / = — A/_i with: 



1 £ 



i i 

= _ 7 + t;Vo 



(31) 



(32) 



and that Eq. l(2"ï"|) thus gives 0, independently of yo- At 
this level of precision we therefore get no information as 
to whether OW's solution is optimal for all y's. 

The next order is "1/A", i.e. we expand both sides 
of Eq. (|22|l and equate terms proportional to —1/A = 
1/ cosh A. We get 



1 f b 

i?i(a)=cosa / Ri(/3)d0. 

2tt J- b 



(33) 



The a dependence on the right hand side is cos a plus a 
constant, so we easily find 



i?i(a) 



sin b 



cosa — 



and 



2 sin b 



7T + b 

The correction to / is given by 

_L f b 

2 sin b 



Vi = 



(34) 



/o = 



cos aRo + Ri (a)] da 



(35) 



7T 



--yi- 



This means that in Eq. I|21|l we get a zeroth order con- 
tribution of 



E, 



GS.O 



2yi -2/i=2/i 



2 sinó 



(36) 



It is easy to see that this expression is the same as the one 
obtained by OW in Ref. jÏÏTJ and if we do a numerical op- 
timization over b we arrive at the notorious 0.434467 . . . 
for the maximal concurrence. This value is obtained for 
b = bow = 1.351802..., corresponding by Eq. (|2*ol to 
y = yow = 0.398316 . . .. 



A. Recursion formula for higher order corrections 

It is tedious, but essentially not difficult to continue in 
the above fashion and calculate higher order corrections. 
A useful trick is to develop a recursion formula. Let us 
write e = 1/1 Al and define: 



R(a)=^R k (a) e k 



dp(a) 
da 



da 



(37) 



The right hand side does not depend on a and we easily 



dO(a,0) _^ de k (a,/3) f 



dft 



d(3 
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The fc'th order terms of Ea . (|22fl give us: 

R k (a) = 

d Pk (a) J_ r b j^d9 k _ 3 (a,f3) 



da 



'j=o 



d(3 



Rj{0) d(3. (38) 



Using the fact that ^3- = 1 we collect terms containing 
Ri. on the left hand side: 



b r 



8(a-j3) + -!- R k {j3) dp 

Z7T 



-EL 

da 



1 

2^ 



-b j=0 °P 



where we have introduced q k {a) as a shorthand nota- 
tion for the r.h.s. The r.h.s. depends only on the known 
functions dp/da and 86/ d (3, and on Rj for j < k. The 
integral operator acting on R k on the l.h.s. of Eq.|j3ÜJ) 
can easily be inverted since it is built from the identity 
and a projection operator (onto a constant). We finally 
end up with the recursion formula 



R k (a) = q k {a) 



1 1 



26 



q k {0)d/3. (40) 



In terms of R k and the auxiliary function q k , we have for 
y k , k > 0: 



y k = 2 [R k (a) - q k {a)} 



(41) 



Note that despite the appearence of a on the right-hand 
side, this relation does make sense since the forní of 
Ea. (|40f) ensures that only terms independent of a sur- 
vive. 

Using Ea. (|4üfl . it is fairly easy to show that 



R 2 {a) 



■ cos 2 a 



sin 6 



cos a 



1 sinó 



26 



cos 6 — 



sin 6 



(42) 



and thus 



2/2 



sin 6 
b + ir 



cos 6 



sin 6 

b + TT 



(43) 



Calculating the first order contribution to the ground 
state energy we find the expression: 



^ 1 b 1 sin6 .,, 

E G s.i(b) = - — - [(6 + 2tt) cos 6 - 2 sin 6 

2 7T 7T + 7T 



(44) 
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Figure 3: Derivative of _Egs w.r.t. e = 1/|A| for e = and 
fixed j/. A positive value for a given 6 indicates that OW's 
solution is a local minimum for the corresponding y. Negative 
vàlues indicates that OW's solution is not a local, and thus 
also not a global minimum. 



In Ea . (|44|l . Bgs.i is expressed as a function of 6, so in 
order to calculate the derivative at fixed y we need to use 
the appropriate implicit differentiation rule. Calculating 
the lowest non-vanishing order we find: 



(dE GS \ 
V de ) 



dE GS t DE y 



Oc 



--E G s,i{b)- 



db 
dE GS0 



v 

2/1 



db dyo j db 
16 6 sin 6 cos 6 
2 7r tt(tt + 6) 



+ 0(e) (45) 



O(e). 



Again we end up with a somewhat complicated expres- 
sion, so we plot its graph in Fig. We note that 
(dE G s/de) y is positive for low 6, but already at 6 = n/2 
(corresponding to y — 1/3) it changes sign and becomcs 
negative. This means that for higher 6's, i.e. lower y's, 
OW's solution cannot be optimal as it is not even a local 
minimum. 

We conclude that in the region of sufficiently large 
magnetizations, i.e. y > 1/3, the OW states (with no 
NN pairs of spins "down") maximize the NN entangle- 
ment locally, i.e. we cannot increase the NN entangle- 
ment by allowing small admixtures of states with NN 
pairs of spins "down". For smaller magnetizations, i.e. 
< y < 1/3, the states that maximize the NN entangle- 
ment necessarily contain NN pairs of spins "down" . 



B. Derivative at fixed y 



C. Higher orders 



As mentioned above, -Egs,i gives us access to whether 
OW's solution is at least a local minimum for a given y. 



The recursion formula l|4U|) is also well suited for nu- 
merical calculations. In Fig. ^ we show a contour-plot 
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Figure 4: í/woif (e, y) calculated numerically from the recur- 
sion formula 14UH . For each b we find the expansion coefficients 
of -ffwoif and y up to the fourteenth order in e = 1/|A|. The 
overall optimal state has y = 0.398316 .... From Sec. IVH Bl 
we know that for all y > 1/3, e = is a local minimum for 
flwoif and this plot indicates that it is also a global one. 



based on such a calculation including all terms up to 
fourteenth order in e — 1/|A|. The plot indicates that 
the calculation in Sec. lVIIBl gives the global answer, i.e., 
for all y > 1/3 the optimal state has no neighboring spins 
"down" . 

Since we perform here the perturbative calculation up 
to the 14th order, we expect that this calculation allows 
us also to obtain some information about the region of 
y < 1/3. From the Fig. 0] (or more precisely from the 
numerical data), one can read off the optimal value of 
e, i.e. optimal value of A. Solving the Bethe ansatz 
integral equation for this value of A we can recover in 
this way the full information about the corresponding 
optimal quantum state. 



VIII. CONCLUSIONS 

In this paper we have studied the question posed by 
O'Connor and Wootters concerning translationally in- 
variant states of N qubits with maximal nearest neigh- 
bor (NN) concurrence. We have answered this question 
for N — > oo using the mapping of the problem onto the 
search for ground states of a certain family of "parent" 
Hamiltonians, described by the XXZ model. Using the 
analytic Bethe ansatz solutions of the XXZ model in the 
limit N — > oo (combining analytic results of low order 
perturbation theory and a numerical calculation of the 
14th order perturbation theory) we have proved that: (i) 



for a given number of spins "down" , i.e. a given mag- 
netization y larger than 1/3, the states that maximize 
the NN concurrence coincide with the ones obtained by 
O'Connor and Wootters, i.e. do not have NN pairs of 
spins "down"; (ii) For small magnetizations, more ex- 
plicitly for < y < 1/3, the states that maximize the 
NN concurrence do contain nearest neighbor pairs of 
spins "down" ; (iii) in particular, the state that maxi- 
mizes the NN concurrence without constraint on y be- 
longs to the family introduced by O'Connor and Woot- 
ters. Our results shed more light on the subtle rela- 
tions between entanglement in spin 1/2 models and the 
ferromagnetic/anti-ferromagnetic character of spin-spin 
interactions. In the appendix we present some simple 
bounds on the optimal magnètic field that corresponds 
to the maximal NN concurrence. 

We acknowledge the support of The Danish Natural 
Science Research Council, DFG (SFB 407, SPP 1078, 
SPP 1116), ESF Program "QUDEDIS", EU FET IST 6th 
Framework Integrated Project "SCALA" , and Spanish 
MEC Grant FIS2005-04627. 



Appendix A: BOUND ON THE OPTIMAL 
MAGNÈTIC FIELD 

If one keeps y fixed, and use Ti. to parameterize the 
s-curve instead of s, one gets: 



dU 



= 2 



dA 1 
~~dM2 ^ 
H 



Vn 2 



i 



-,p 



n 



Vh 2 



i 



-i y 



where P is the probability of two neighboring spins being 
both "down" . 

Demanding that d ^jf = 0, we find: 

H opt = ±y/(2VPy/P + y) 



where 7í pt is the optimal magnètic field. 
Hopt = iff y = Q. 

Using a simple bound on P, 



Obviously, 



P < (y + l)/2 in the limit N -> oo, 
we obtain a lower bound on 7í op t'. 



> 2//(Vy + iV3y + i). 

This bound does not work well for y > 1/3, because it 
gives a finite bound, maximized for y = 1 when we find 
Ti-opt > l/(2\/2), whereas we know that in this regime of 
y's 7í pt = 00 ■ For smaller vàlues of y < 1/3, both the 
optimal e (i.e. A, see Fig.0J, as well as the optimal 7i pt 
attain finite vàlues, so that the bound might become more 
useful. In particular, the results of Fig. 0] suggest that 
as y approaches zero, the optimal e approaches 1 more or 
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less linear ly, as 1 — (1/3) 1 y, which in turn implies that 
the optimal A approaches — 1 as — 1 — 3y. Thus for small 
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